Neel to valence-bond solid transition on the honeycomb lattice: Evidence for 

deconfined criticality 
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We study a spin-1/2 Heisenberg model on the honeycomb lattice with nearest-neighbor anti- 
ferromagnetic exchange J that favors Neel order, and competing 6-spin interactions Q which fa- 
vor a valence bond solid state in which the bond-energies order at the "columnar" wavevector 
K = (27r/3, — 27r/3). Using projector Quantum Monte Carlo techniques, we find evidence for a 
direct continuous quantum phase transition between the Neel and VBS states, analogous to the 
deconfined critical point between these states on the square lattice. This implies that such decon- 
fined critical points can exist on the honeycomb lattice although Berry phase effects allow tripled 
hedgehog defects of the Neel order parameterto play a role in the critical theory, unlike the square 
lattice case in which only quadrupled hedgehog defects are allowed. 

PACS numbers: 



Many interesting materials appear at low temperature 
to be on the verge of a quantum phase transition in- 
volving a qualitative change in the nature of the ground 
state[T]. When one of the two competing T — Q phases 
spontaneously breaks a symmetry, the transition can 
be studied using a path integral representation with a 
Landau- Ginzburg actional written in terms of the or- 
der parameter that characterizes the broken symmetry 
phase[T]. If phases on two sides of the critical point break 
different symmetries, Landau- Ginzburg theory generi- 
cally predicts a direct first-order transition or a two- 
step transition with an intermediate phase. However, 
this path integral description in terms of order-parameter 
variables can sometimes involve Berry phases in a non- 
trivial way[3H5]. The presence of Berry phases, which 
correspond to complex Boltzmann weights for the corre- 
sponding classical statistical mechanics problem in one 
higher dimension [T], can invalidate some of the conclu- 
sions reached by the Landau-Ginzburg approach. 

In some of these cases, it is useful[6] to think in terms 
of topological defects in one of the ordered states, and 
view the competing ordered state as being the result 
of the condensation of these topological defects — this 
description [6j makes sense only if the quantum numbers 
carried by defects in one phase match those of the or- 
der parameter variable in the other phase. Under cer- 
tain conditions, this alternate "non-Landau" description 
generically predicts a direct continuous transition^ [S] 
between the two ordered states, in contrast to predic- 
tions of classical Landau-Ginzburg theory. Square-lattice 
quantum antiferromagnets undergoing a transition from 
a ground state with non-zero Neel order parameter to 
a valence-bond solid (VBS) ordered state, in which the 
"bond-energies" (singlet projectors) P^^) = \ — Si ■ Sj 
on nearest-neighbour bonds {ij) in the x (y) direction 
develop long-range order at the "columnar" wavevectors 
Ki = (tt, 0) (K2 — (0, tt)), provide the best-studied ex- 
ample of such a non-Landau transition. In this case, Z4 



vortices in the VBS order carry a net spin S = 1/2 in 
their core, suggesting that the onset of Neel order can be 
studied using a CP^ description of A/^: Mg — z'^UafiZp, 
where a are the three Pauli matrices and the spinful Z4 
vortices are represented by a two component complex 
bosonic field coupled to a compact U(l) gauge field 
yi^j 6 . Since Berry phases associated with space-time 
monopoles of ^;i[3Hil IS] can be argued[71 [S] to lead to 
a complete suppression of these monopoles at criticality, 
the Neel- VBS transition is described within this picture 
by a 2+1 dimensional non- compact CP^ (NCCP^) theory. 
Although some violations [TDHIl] of conventional scaling 
behavior need to be better understood, this prediction of 
a "deconfined critical point" for the Neel- VBS transition 
on the square lattice is largely consistent [TTl fTSHlQ] (see 
also Ref. [20]) with numerical data at sizes accessible to 
Quantum Monte Carlo (QMC) simulations of computa- 
tionally tractable spin models. 

Here, we study a spin-1/2 Heisenberg model on the 
honeycomb lattice with nearest-neighbor antiferromag- 
netic exchange J that favors Neel order, and compet- 
ing 6-spin interactions Q which favor VBS order at the 
columnar wavevector K = (27r/3, — 27r/3): 
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where {{ijklmn)) denotes hexagonal plaquettes (Fig[T]). 
We use a projector QMC method pTll^ to study the cor- 
responding r = Neel- VBS transition, and obtain evi- 
dence for a direct continuous transition at {Q/ J)c = Qc ~ 
1.190(6), with correlation length exponent ly w 0.54(5), 
and anomalous exponents fyNcci ~ 0.30(5), and jyvBS ~ 
0.28(8). We suggest that this critical point also describes 
universal properties of the transition to plaquette VBS 
order recently observed in the J1-J2 model on the 
same lattice. 

The continuous nature of the transition studied here 
is in sharp contrast to the strongly first order transitions 
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to staggered VBS states observed on square and honey- 
comb lattices [261 [27]. To understand the significance of 
our result, we note that Z-^ vortices in this VBS state on 
the honeycomb lattice carry a net spin S" = 1/2 in their 
cores (Fig[l]), pointing again to a CP^ description of the 
corresponding transition. In this description, space-time 
monopoles in A^j, correspond to hedgehog defects of the 
Neel order parameter (H j^. The creation operator for 
these monopoles transforms under lattice symmetries in 
the same way as the order parameter ^ for VBS order 
at wavevector K = (27r/3, —2tt/3) on the honeycomb lat- 
tice and wavevectors K1/2 on the square lattice, allowing 
one to view these VBS states as condensates of these 
hedgehogs [H [3 [a H]. 

In particular, the monopole creation operator picks up 
a 27r/3 phase under honeycomb lattice rotations. There- 
fore, insertions of tripled monopoles are allowed in the 
CP^ description of the honeycomb lattice transition to 
VBS ordering at wavevector K. If these are relevant, one 
expects the correct long-wavelength description to be a 
conventional Landau-Ginzburg theory written in terms of 
the Neel order parameter Ms and the VBS order param- 
eter ^, and the transition to be first-order in the simplest 
scenario [3 [5], or proceed in two steps with an interme- 
diate phase. Thus, the direct second-order transition ob- 
served here strongly suggests that tripled monopoles are 
irrelevant perturbations at the T = critical point of 
the two-dimensional NCCP^^^ theory for N = 2. This 
result is difficult to anticipate in the absence of numeri- 
cal evidence, since the physical N — 2 case lies between 
two contrasting extremes: At = 1, tripled monopoles 
are relevant [SJ [71 1281 129) and lead to a weakly-first or- 
der transition\30\ . while in the N 00 limit, they are 
strongly irrelevant [5l 171 [28] . In contrast, in the original 
square lattice case, 4* picks up a phase of n/2 under lat- 
tice rotations, and only quadrupled monopoles are allowed 
in the CP^ description. These can be argued]?] [S] to be 
irrelevant at the physical N = 2 critical point by noting 
that they are irrelevant both at TV = 1 [5] [7] [28 ] [29 ] l3T]. 
and in the iV — ^ 00 limit [5] [7] [5S], suggesting that they 
are also irrelevant in the physical N — 2 case. 

All of the foregoing suggests that the transition from 
the Neel state to a plaquette VBS state ordered at 
wavevector K (Fig[T]), first observed in exact diagonal- 
ization [23] of a frustrated Ji — J2 spin 1/2 model on 
the honeycomb lattice, is also a second-order transition 
in the universality class of NCCP^ critical point. Indeed, 
some indications of a continuous transition have recently 
been obtained using DMRG methods [24] [25]. 

We study H on L x L honeycomb lattices (Fig [T]) of 
2L^ spins, with periodic boundary conditions and L a 
multiple of 12 up to L = 72. We use a T = projector 
QMC algorithm [22], with a sufficiently large projection 
length m = cL^ (with c ranging from 4 to 12) to ensure 
convergence to the ground state — all expectation val- 
ues (...) in what follows thus refer to the singlet ground 




FIG. 1: (Color online) The honeycomb lattice has a two-site 
basis (labeled A and B) and elementary Bravais lattice trans- 
lations ei and 62, with distances from origin specified in units 
of ei and 62. Three types of bonds (labeled 0, 1, 2), oriented 
along the three principal directions, "belong" to each Bra- 
vais lattice site. If dimers with solid (hashed) fill pattern are 
taken to denote high (low) values of {P(ij)) on the correspond- 
ing link (ij), the sample is in a columnar VBS state, ordered 
at wavevector K, with phase of the order parameter pinned 
to 47r/3, 0, and 2n/3 in black, brown and blue regions respec- 
tively. Conversely, if dimers with solid (hashed) fill pattern 
are taken to denote low (high) values of {P{ij)}, the sample 
is in a plaquette VBS state, ordered at wavevector K, with 
the phase of ^ pinned to tt/3, tt, and 5n/3 in black, brown 
and blue regions respectively. Three domain walls meet at 
the core of the Z3 vortex, which carries a net spin of 5 = 1/2. 
Also shown is a pictorial depiction of the six-spin interaction 
terms that drive the transition to VBS order. 

state of H on the L x L lattice. At small q, the ground- 
state is Neel ordered, as characterized by the Neel order 
parameter Mg = ^ where m is the local Neel 

order parameter field: 

rn{r) = S^a ~ SpB ■ 

Here fA (fB) refers to the A {B) sublattice site belonging 
to Bravais lattice site r (Fig. Un and the S are spin-1/2 
operators. To characterize the strength of Neel order, we 
compute (Mg) and (M^), as well as correlations ofm{r), 

C{RaL) - ^ ^{m{r) ■ m{r + R^l)) , 

r 

at separations RaL = {aL,aL) with a = ^, ^, ^, |. 

To locate the transition at which Neel order is lost, 
we use two "dimensionless" quantities. The first is the 
Binder cumulant g^g = {{M^)'^) / {M'^y , while the sec- 
ond is a ratio of Neel correlators at different distances, 
for instance C{Rl/^)/C{Rl/q). Both quantities are ex- 
pected to obey a scaling form Fg ^c[^<1n) if Neel or- 
der is lost via a continuous transition at qcN- Here, 
Fg^j ,c &re universal scaling functions of the argument 
AqN = {q — qcN)L^^^" and vn is the correlation length 
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FIG. 2: (Color online) Binder cumulant of Ms as a func- 
tion of q for various L (symbols), fit to a polynomial in 
(g — qcN)L^^"'^ (lines) using best-fit values vm = 0.5080 and 
QcN = 1.1912. Inset : q dependence of the Neel correla- 
tion ratio C(i?L/3)/C(i?i/6) for various L (symbols), fit to 
a polynomial in [q — qcN)L^^'^" (lines) using best-fit values 
UN = 0.4984 and qcN = 1.1944. All best-fit values are for the 
L > 48 part of the displayed data. 



FIG. 3: (Color online) Binder cumulant of as a function 
of q, for different sizes L (symbols), fit to a polynomial in [q — 
qcD)-L^^'''^ (lines) with best-fit values ud = 0.5237 and qcD = 
1.1892. Inset : VBS correlation ratio D{Rl/3) / D{Rl/6) for 
different L as a function of q (symbols), fit to a polynomial in 
{q — qcD)-L^^'''^ (lines) with best-fit values vd ~ 0.5721 and 
qcD = 1.1844 (fines). All best-fit values are for the L > 48 
part of the displayed data. 



exponent associated with Neel correlations. Both quanti- 
ties allow us to locate a continuous Neel ordering transi- 
tion by looking for a crossing point of g-dependent curves 
for various L. In the vicinity of such a transition, we 
expect the scaling forms (aI^) = L-(i+''N<=<=i)G'jg (A^at), 

C{RaL) = L-(i+''~-i)Gc(AgAr) for the corresponding di- 
mensionful quantities. 

At large g, we find that VBS order develops at the 
columnar wavevector K. We characterize this in terms 
of the VBS order parameter ^ — ^ ^ - Vp, where Vf: is 
the local VBS order parameter field: 

Here, Pf^ (^ = 0,1,2) denotes the singlet projector on 
the bond n "belonging" to Bravais lattice site r (Fig.[T]). 
As emphasized in Fig. [TJ there are two distinct possi- 
bilities for three-fold symmetry breaking VBS order at 
wavevector K. These columnar and plaquette VBS or- 
derings (Fig [T]) can be distinguished by the phase of 
the VBS order parameter ^. In QMC simulations, in- 
formation about the phase of ^ can be obtained from 
the estimator E^, whose average E^^, over the QMC run 
gives the quantum-mechanical expectation value {^) . Al- 
though E^j, is a basis dependent quantity, the histogram 
of its phase can nevertheless be used to distinguish pla- 
quette VBS states from columnar states at the same 
wavevector|15[ [T7]. To quantify the strength of VBS or- 
der without distinguishing between these two possibili- 
ties, we compute (5'^^^') and correlations of Vp, 



at the same separations RaL as for Neel correlators. 

To locate the VBS transition, we again focus on 
two dimensionless quantities. The first is the (basis- 
dependent) Binder cumulant[32 of i?^ defined as gE^, = 

|£'^|'*/ , and the second is the ratio of the VBS 

correlators at different distances, say D{R]^/^)/D{R]^/q). 
These quantities are again expected to obey scaling forms 
Fgj^^ D{AqD) if VBS order is lost via a continuous tran- 
sition at QcD- Here, Fg^^o are universal scaling func- 
tions of the argument Ag^ = {q — qcD)L^^'^'^ and 
is the correlation length exponent associated with VBS 
correlations. A crossing point for different L as a func- 
tion of q for each of these two quantities then reveals 
a continuous VBS ordering transition. Close to such a 
continuous transition, we also expect the corresponding 
scaling forms (l^-p) = i-(i+''™s)G^(Agi5), D{RaL) = 
^-(i+i)vBs)(^^^^g^'j fQj. dimensionful observables. 

Our analysis uses the crossings of the Binder ratios 
and gEq, to locate the Neel and VBS transitions — at this 
stage, we do not assume that the two transitions coin- 
cide. Given the relatively sharp nature of the crossings 
(displayed in Fig [2] and Fig[3]) and the monotonic nature 
of their q dependence for fixed L, we are confident that 
the transition(s) is (are) continuous. The behaviors of 
the Neel and VBS correlation ratios in the critical region 
are also shown as insets, and confirm this finding. 

Therefore, we fit the data for each dimension- 
less quantity (g^g , and the correlation ratios 

C{RL/3)/CiRL/6) and D{Rl/s)ID{Rl/^)) in the criti- 
cal range to a polynomial function of {q — qc)L^^'^ (cor- 
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FIG. 4: (Color online) Scaling collapse of Binder ratios of 
Ma and , and of correlation ratios C{Rl/3)/C{Rl/q) and 
D{Rl/3)/D{Rl/q), using values of q^n, QcN, vd and 
quoted in legends of Fig. [2] and Fig. [3] 
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FIG. 5: (Color online) Scaling collapses of (Mf), 
C(L/3,L/3) and VBS D(L/3,L/3) , obtained using the fol- 
lowing best-fit values: Qcn ~ 1.1956, = 0.5003, r^Neei = 
0.3539 {{M!)),qcN = 1.1944, un = 0.5085, 77N001 = 0.3066 
(C(i?i/3)), gcD = 1.1864, i/D = 0.558, r?vBS = 0.25 ((|*|'')), 
and q^D = 1.1848, i^d = 0.5733, r;vBS = 0.2961 {D{RLfi)). 
Quoted best-fit values are for the L > 48 part of the displayed 
data. 



responding to a polynomial approximation of the corre- 
sponding scaling function), with the corresponding q^, 
V and the polynomial coefficients being fitting parame- 
ters. For each dimensionless quantity, the best-fit values 
vary somewhat depending on the range of L and q stud- 
ied; results of such a fit for one choice of data-set are 
displayed as the lines in Fig. [2] and Fig. |3] The cor- 
responding scaling collapse of each dimensionless quan- 
tity is displayed in Fig. [4j In addition, we fit each of 



the dimensionful quantities (M^), C{Rl/3), (I^'P), and 
D{Rl/3) in the critical region to times a poly- 

nomial in {q — qc)L^^'^ ■ Apart from qc, v and polynomial 
coefficients, these fits involve the appropriate anomalous 
exponents r/Ncci or ?7vbs as additional fitting parameters. 
For each dimensionful quantity, the best-fit values again 
vary somewhat depending on the critical range (in L and 
q) used for the fit. Fig [s] displays the scaling collapses 
obtained from best-fit values for one choice of data-set. 

Based on a detailed study of such fits, we estimate 
q,N ~ 1.1936(24), q,D « 1.1864(28), m = 0.51(3), 
lyjj = 0.55(4), ?7Noci = 0.30(5) and ?7vbs = 0.28(8). The 
error bars quoted here reflect not just the error in deter- 
mining best-fit values for a given data-set for each quan- 
tity, and variation in these best-fit values from quantity 
to quantity, but also the dependence of these best-fit val- 
ues on the data set used, i.e. the size of the critical 
window in q, and the range of L used in the fits. We also 
emphasize that our estimates of ?7vbs and fyNoci depend 
sensitively on the value of Qc, resulting in the relatively 
large error bars quoted here. Nevertheless, we are in a 
position to exclude the relatively tiny values of 77 that 
characterize conventional second-order critical points in 
2-1-1 dimensions. Since i^^v coincides with Pjj within er- 
ror bars, and the allowed ranges of Qcn and qco almost 
touch at the one-sigma level, the simplest interpretation 
of our data is that Neel order is lost and VBS order sets 
in at a single continuous transition whose location is es- 
timated to be Qc ~ 1.190(6), with correlation exponent 
v = 0.54(5), and anomalous exponents TyNcoi — 0.30(5) 
and 77VBS — 0.28(8). This, taken together with the 
relatively large values of ryNooi and t^vbs characteristic 
of deconfined critical points, suggests an interpretation 
in terms of deconfined criticality. Indeed, our estimate 
of the universal critical value g* = 1.42(1) of the Neel 
Binder ratio at is consistent with the value for the 
same quantity quoted in Ref [T7| for the square lattice 
transition. Further, our estimates of ryvBS and 77nco1 are 
compatible within our errors with the best available esti- 
mates from square lattice studies P^. However, we note 
that our estimate of i' is slightly lower than the best pre- 
vious estimate of ~ 0.59 at the Neel- VBS transition on 
the square lattice [TT]. 

In summary, we find a direct continuous transition be- 
tween a Neel state and a VBS state ordered at the colum- 
nar wavevector K on the honeycomb lattice, which ap- 
pears to share universal properties with deconfined phase 
transitions on the square lattice. We suggest that this 
critical point also describes universal properties of the 
transition to plaquette VBS order recently observed [531- 
[25] in the J1-J2 model on the same lattice. Our results 
motivate follow-up work that explores [33] possible loga- 
rithmic violations of standard scaling for some quantities; 
this may elucidate the nature and likely origins of such 
violations seen earlier on the square lattice [TU1 - IT2] . 
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